Quantum Monte Carlo study of large spin-polarized tritium clusters 
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This work expands recent investigations in the field of spin-polarized tritium (TJ,) clusters. We 
report the results for the ground state energy and structural properties of large T4. clusters consisting 
of up to 320 atoms. All calculations have been performed with variational and diffusion Monte 
Carlo methods, using an accurate ab initio interatomic potential. Our results for A'' < 40 are 
in good agreement with results obtained by other groups. Using a liquid-drop expression for the 
energy per particle, we estimate the liquid equilibrium density, which is in good agreement with our 
recently obtained results for bulk TJ,. In addition, the calculations of the energy for large clusters 
have allowed for an estimation of the surface tension. From the mean-square radius of the drop, 
determined using unbiased estimators, we determine the dependence of the radii on the size of the 
cluster and extract the unit radius of the T4, liquid. 

PACS numbers: 67.65.-|-z,02.70.Ss 



o 



> 

(N 
(N 



O 
O 



X 



INTRODUCTION 

The extreme quantum nature of electron spin- 
polarized hydrogen (HJ,) and its isotopes, spin-polarized 
deuterium (DJ,) and spin-polarized tritium (TJ,), pro- 
moted renewed theoretical interest [IHll in their con- 
densed phases. They are characterized by the small mass 
of their atoms and their weakly attractive interatomic po- 
tential, which is very accurately determined in ah initio 
calculations. Already in the seventies, several theoret- 
ical predictions of the spin-polarized hydrogen systems 
appeared 0,[1]- In 1976, Stwally and Nosanow under- 
lined HJ, as a first candidate for achieving a Bose-Einstein 
condensate (BEC) state. Their theoretical prediction was 
experimentally confirmed in 1998 by Fried et al. [lo| . who 
succeeded to overcome demanding experimental obsta- 
cles and reported the formation of a BEC state with HJ, 
atoms. 

Further investigation related to possible candidates for 
BEC state in hydrogen systems was done by Blume et al. 
In their work [ll| , spin-polarized tritium clusters as well as 
optically pumped tritium condensate were theoretically 
investigated for the first time. Blume et al. reported their 
DMC results for the ground state energy and structural 
properties of TJ, clusters consisting of up to 40 atoms. In 
addition, it was also shown that the smallest TJ, cluster 
is a trimer, i.e. (T4,)3; negative ground-state energy was 
not obtained for the dimer (T4,)2, which means that TJ, 
trimer is an example of Borromean or halo state. The 
same conclusion for (T4,)3, obtained with the finite ele- 



ment method, was reported by Salci et al ll|. Because of 
evident resemblance of bosonic TJ, and ^He atoms Blume 
et al. [H also compared general properties of both types 
of clusters. They showed that common attributes of TJ, 
clusters are weaker binding and greater interparticle dis- 
tances between atoms, in comparison with "^He clusters 
having the same number of atoms. In the same work, 
results of coupled-channel scattering calculations for two 
TJ, atoms are reported, indicating the possibility for for- 
mation of a tritium condensate using its broad Feshbach 
resonance. 

Mixed clusters consisting of spin-polarized hydrogen- 
tritium and deuterium-tritium atoms have also been in- 
vestigated. 0, m It has been shown that three T\. atoms 
are needed to bind one DJ, atom in a stable system, mak- 
ing thus clusters (T4,)7vDi stable for all > 3. On 
the contrary, it has been shown that even 60 TJ, are not 
enough to bind one HJ, in a stable system. Namely, the 
ground-state energy of the cluster (T4,)6oH4, is within the 
error bar equal to the ground state energy of (T4,)go, lead- 
ing to the conclusion that clusters (T4,)jvH4, for A^ < 60 
are effectively unstable or are at the threshold of binding. 
Due to the more complicated calculations in the case of 
several fermionic DJ, atoms, so far only the stability limits 
of small mixed spin-polarized deuterium-tritium clusters 
having up to 5 DJ, atoms have been examined. P] 

Despite the lack of experimental verification, bulk 
properties of all spin-polarized hydrogen isotopes have 
been theoretically predicted. Bulk properties of the D\. 
system are conditioned by the number of occupied nu- 
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clear spin states. [l2-14| If only one nuclear spin state 
is occupied (D4,i), the system is in a gas state at zero 
pressure, while in the case of two (0^,2) or three (04,3) 
equally occupied nuclear spin states, the system remains 
liquid at zero pressure and zero temperature. Exten- 
sive investigations of the HJ, and TJ, bulk systems have 
been carried out recently with the diffusion Monte Carlo 
(DMC) method 0,11], which provides exact results within 
errorbars for bosonic systems. Using the DMC method 
for bulk K-l and T-l, the energy per particle, structural 
properties, as well as densities and the pressure of the 
gas (liquid)-solid transition have been predicted. An ac- 
curate calculation of the ground-state energy per particle 
in bulk HJ, has allowed the confirmation of its gas nature 
in the limit of zero temperature and up to 170 bar, point 
at which HJ, solidifies. A similar investigation of bulk 
T], has revealed that the system is a liquid up to 9 bars, 
where it crystallizes. 



In this work, we expand previously reported studies of 
pure TJ, clusters. [H We report the ground-state en- 
ergy of clusters having up to 320 atoms, as well as their 
structural properties, obtained with the DMC method. 
From the density profiles, we estimate the thickness of 
the clusters' surface. Justification for carrying out de- 
manding calculations for large clusters lies in the fact 
that the present results for clusters can be used to extrap- 
olate precise equilibrium TJ, bulk properties. A goal of 
our investigation is to examine the validity of the liquid- 
drop formulas when they are applied to TJ, clusters, as 
it was done in the past for ^He and ^He clusters. [23.l25j 
In helium clusters, liquid-drop formulas were successfully 
applied and the results for the equilibrium energy per 
particle and the unit liquid radius were in good agree- 
ment with experimental studies. [13, [111 We have used 
the energy per particle of the TJ, clusters to extrapo- 
late the equilibrium energy per particle in bulk TJ,. We 
compare the estimated result with the energy perparticle 
calculated in a recent DMC study of a bulk T^. l&i In ad- 
dition, the surface tension of liquid TJ, is estimated and 
compared with known results for ^He and ^He liquids. 
Furthermore, we extract the unit radius of the liquid us- 
ing the average distance of the particles to the centre of 
mass of the cluster. 



In Sec. II, we report briefly the DMC method and 
discuss the trial wave functions used for importance sam- 
pling of the clusters. Sec. Ill reports the results obtained 
by the DMC simulations. Finally, Sec. IV comprises a 
summary of the work and an account of the main con- 
clusions. 



METHOD 

The starting point of the DMC method is the 
Schrodinger equation written in imaginary time, 

-h^^^ = iH-E,)^iR,t) , (1) 

where E-^ is a constant acting as a reference energy and 
R = (ri, . . . , rjv) collectively denotes particle positions. 
The iV-particle Hamiltonian, H, is given as 
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where V{r) is the interaction potential. The inter- 
atomic interaction between tritium atoms is described 
with the spin- independent central triplet pair potential 
j„ , which was determined in an essentially exact way 



by Kolos and Wolniewicz [15|. As in our recent DMC 
calculations of bulk HJ. and T4,,0, [1] we have used the 
recent extention of Kolos and Wolniewicz data to larger 
intcrparticle distances by Jamieson et al. (JDW). [16[ 
The potential is finally constructed using a cubic spline 
interpolation of JDW data, which is smoothly connected 
to the long-range behavior of the T4,-T4, potential as cal- 
culated by Yan et al. The JDW potential used in the 
present work has a core diameter a — 3.67 A and a mini- 
mum of —6.49 K at a distance 4.14 A. We have previously 
verified that the addition of mass-dependent adiabatic 
corrections (as calculated by Kolos and Rychlewski [l8| |) 
to the JDW potential does not change the energy of bulk 
spin-polarized tritium. Q It is worth mentioning that 
within the Born-Oppenheimer approximation it has been 
explicitly shown that in the spin- aligned electronic state, 
tritium nuclei behave as effective bosons, [l^ 

DMC solves stochastically the Schrodinger equation 
(HJ) by multiplying with the tp{R), a trial wave 

function used for importance sampling, and rewriting 
Eq. (HI) in terms of the mixed distribution i) — 

t)V'(-R). Within the Monte Carlo framework, 
^{R,t) is represented by a set of walkers. In the limit 
t — > 00 only the lowest energy eigenfunction, not or- 
thogonal to ip{R), survives and then the sampling of the 
ground state is effectively achieved. Apart from statisti- 
cal uncertainties, the energy of a A^-body bosonic system 
is exactly calculated. 

In the present simulations Jastrow trial wave functions 
have been used. 
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with a two-body correlation function /(r) 
f{r,j) = exp 



(3) 
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where n = 1 or n = 2, depending on the size of the 
cluster, and b and s are variational parameters. Previ- 
ous experience in work with small pure and mixed spin- 
polarized tritium clusters Q has shown that the best 
choice for the two-body correlation function is obtained 
with n — 1. We have used this type of function for clus- 
ters having < 60 atoms, but for larger clusters the 
variational energies obtained with this type of function 
worsen significantly. Then, for larger clusters, it is better 
to consider the model with n — 2, which is similar to f{r) 
which has been used in recent investigations of vortices in 
large ''He clusters. [2^ Eq. ^ with n — 2 defines much 
better the confinement of space in which large number of 
TJ, atoms is settled and definitely provides better VMC 
energies for clusters having > 80 atoms. 

The optimization of the trial wave functions has been 
done for all clusters by means of the variational Monte 
Carlo method. For clusters with iV < 60 the best varia- 
tional parameters vary from 6=3.574 A to 6=3.605 A and 
from s=0.0328 A"! to s=0.0073 A"! for increasing N. 
In case of clusters having > 80 atoms, the parame- 
ter b assumes values from 3.574 A to 3.605 A, while at 
the same time s varies from 000162 A^^ 0.0000145 
A^^. It is worth noticing that in both types of two-body 
correlation functions (n = 1,2) the parameter 6 remains 
practically constant, while s always decreases with A^. 

For several clusters we have verified that 1000 walkers 
are enough for excluding the bias coming from the size 
of the population ensemble used in a simulation. Thus, 
we have decided to employ this number of walkers in all 
the remaining DMC calculations. The same conclusion 
emerged also from previous experience in pure and mixed 
T| clusters. 

In order to eliminate bias coming from the time-step 
value used in simulations, all calculations have been per- 
formed with several At time-steps which assume values 
within the interval 5 x 10""- 1.3 x lO'^K'^ From the ob- 
tained results, we have extrapolated the result to At 0. 
In accordance with the DMC method used in this work, 



which is accurate to second order in the time step [2l| , the 
extrapolation is made with a quadratic function. Second 
order DMC enables the use of greater time-steps than 
the linear DMC method. 



RESULTS 

Our DMC results for the ground-state energy per par- 
ticle and radii of the investigated clusters are given in 
Table HI For clusters consisting of up to 10 T| atoms 
we have already shown j^lgood agreement with results 
obtained by Blume et al. Ij. Here, we extend this com- 
parison for clusters consisting of up to 40 T| atoms. In 
Fig. [U comparison of our results and the ones by Blume 
et al. is shown. As in the case of small clusters, we 
report slightly lower ground-state energy for all clusters 



(TJ,)^?. These differences are mainly due to the fact that 
there are small differences in the potential of interaction 
employed in two simulations. Namely, Blume et al. [ij 
have included in the Hamiltonian the damped three-body 
Axilrod- Teller potential term [22] which causes a slight 
raise of the ground-state energy. 

In Ref. T, weaker binding and greater spread in the TJ, 
clusters were emphasized as the main difference between 
small (TJ.) AT and ^He) n clusters, for the same A^ and up 
to 40. Similar comparison can be done for large clusters. 
Namely, the energy per particle of the largest investi- 
gated "^He cluster with the DMC method, (''He)ii2, is 
-3.780(3) K [2^. This can be compared with the energy 
per particle of the largest TJ, cluster, (T4,)32o, which is - 
2.286(8) K (Table[l|. It is clear from these resuhs that in 
the ("He)ii2 cluster, which consists of almost three times 
smaller number of atoms than (T4,)32o cluster, binding is 
significantly stronger. From that, we conclude that the 
binding in large TJ, clusters is very weak, as it was al- 
ready concluded for small clusters [l|. 

The energy per particle of quantum liquid clusters 
as a function of A^ is well reproduced by a liquid-drop 
model, (Uliil 



E{N)/N = Ey + xEs 



x^E, 



(5) 



where Ey, Eg and Ec are respectively the volume, surface 

and curvature terms, and the variable x is defined as 
x=A^-i/3. 

In Fig. [21 we have plotted results for the energy per 
particle from Table [H as well as a line on the top of 
data, which represents the best fit with expression 
We have included in our fit all the investigated clusters 
and the best set of parameters obtained with the above 
mentioned fit is: E^^ -3.66(3) K, £;,= 10.2(2) K and 
E,= -6.1(4) K. 

The parameter E^ represents the energy per particle of 
bulk liquid TJ, at the equilibrium density. This extrapo- 
lated result is in a very good agreement with our recent 
results obtained in calculations of bulk TJ, @ . With the 
DMC method we obtained po = 0.007466(7) A'^ as the 
equilibrium density of liquid TJ, and eo = —3.656(4) K 
as the energy per particle at that density. 

It is also important to emphasize that the parameters 
Eg and E^ have been obtained without including in 
the fit the bulk energy per particle at equilibrium. Fur- 
thermore, the parameters in Eq. ([5]) remain practically 
the same when the energy per particle at the equilib- 
rium density is included in the fit. Thus, the demanding 
DMC calculations have been worthwhile because the re- 
sult of our liquid-drop model ([5]) does not depend on 
the knowledge of the equilibrium energy per particle of 
the bulk. Also, we have determined that (TJ,)28o is the 
'smallest' cluster needed to be included in the fit in or- 
der to extrapolate the parameter Ey properly. Namely, 
using just DMC results for clusters having A" = 20 — 280 
atoms in the fit we get Ey= -3.69(3) K, which is within 
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the error bars the same as the result obtained inchiding 
the largest cluster (iV=320, £'„=-3.66(3) K). However, 
the extrapolation of the equilibrium energy per particle 
with fits including just results for clusters smaller than 
(T4,)28o always produces lower energy than the one cal- 
culated for the bulk. For example, the extrapolated bulk 
energy with results for clusters having up to 120 atoms 
is -3.85(2) K, which is around 5% lower than eg. Sim- 
ilar conclusions about this fit emerged for ^He clusters 
[25| . where it was emphasized that an accurate extrapo- 
lation of the bulk equilibrium energy from finite cluster 
calculations should include relatively large clusters. 

We have also tried to fit all the obtained data for E/N 
with a linear function in iV~^/^, as in Ref. l2i, but that 
kind of fit has not been so precise as the one performed 
including a quadratic dependence on N~^^^. With the 
linear fit E^— -3.26(3) K, which is around 11% higher 
than the ep obtained for bulk, showing the necessity of 
including the second-order term. 

The second parameter Eg extracted from Eq. ^ is 
related to the surface tension of liquid TJ, through 



and 
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where tq is the unit radius of the liquid. The unit radius 
of the liquid can be determined in two ways, using the 
result of the equilibrium density of the bulk liquid 
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or from the expression 



ro(iV) = 
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(7) 



(8) 



where (r^ {N)) is the mean-square radius of a cluster with 
N atoms. 

Using the result for the equilibrium density of the TJ, 
liquid po = 0.007466(7) A'^ [g] and Eq. we obtain 
ro =3.18(1) A . 

The second method for obtaining the unit radius, 
Eq. @, was previo usly used in the study of "^He clus- 
ters. (2J,|25| In Ref. |2J, it is emphasized that only those 
^He clusters with more than ten atoms have a radius 
which increases approximately as N^^^ . Since in our cal- 
culations we obtain unbiased mean-square radii of clus- 
ters with pure estimators, we tried to interpolate our 
data for clusters radii, Rcm{N) — y/ {r^{N)), with sev- 
eral polynomial functions of the variable x = N-^^^. As 
in the case of the interpolation of the energy per particle, 
we have included all clusters having = 20 — 320 atoms. 
In Fig. 131 clusters' radii obtained from the calculations 
are plotted, as well as two lines on top of data which 
represent two interpolations, using functions 
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The parameters extracted from the fits are: a = 
1.85(12)A, b = 2.34(3)A, c = 2.55(1)A and d = 
3.73(12)A. Using these interpolation parameters and the 
definition of unit radii given in (|8]), in the limit N ^ oo, 
we extract an equilibrium radius rg = 3.02(4) A using 
the function (jH) and tq = 3.29(1) A using the function 
(fTOl) . Since in both cases the quality of the fit is very 
good, we cannot state which of the two extracted results 
for equilibrium unit radius should be considered as the 
better estimation. We can thus only conclude that the 
equilibrium unit radius assumes a value within the inter- 
val from 3.02(4) A to 3.29(1) A. 

Therefore, the estimation of the equilibrium unit ra- 
dius from the results obtained in clusters calculations 
is very sensitive to the choice of interpolating function. 
Because of that, we decided to include the unit radius 
3.18(1) A , derived from the bulk TJ,, in the estimation 
of the surface tension. 

It is useful to compare the value rg with a = 3.67 A. 
If we consider rg as the radius of the sphere that one TJ, 
atom occupies at equilibrium density in a liquid, a has 
to be smaller than 2 rg, as it is in our case. On the other 
hand, the value of rg of liquid TJ, is also greater than the 
unit radius of hquid ^He (ro=2.1799 A) This also 

explains the greater spread in TJ, clusters because it is 
clear that TJ, atoms occupy more space than ^He atoms 
at equilibrium density. 

From ro, it is possible to calculate the liquid surface 
tension t using expression ([6]); we have obtained t = 0.08 
kA~2. There is no experimental result for the surface 
tension of liquid TJ,, and this prediction is, to the best 
of our knowledge, the first estimation of t for liquid TJ,. 
Contrary to liquid TJ,, the surface tensions of "^He and 
■^He liquids have been experimentally investigated and 
the measured values are respectively 0.27 KA~^ and 0.11 
KA-2. [23 In the case of liquid TJ, our estimated value 



bx. 



(9) 



of the surface tension is even smaller than the surface 
tension of the '^He liquid, although bulk TJ, is a bosonic 
system. Explanation for such a small value of the surface 
tension lies in the fact that the interaction between TJ, 
atoms is described with a very shallow potential. 

In addition to the ground-state energy, we have also 
studied the structure of TJ, clusters. Exact estimators 
of DMC method have been employed to calculate values 
such as the pair distribution function P{r), as well as 
the distribution of particles with respect to the centre 
of mass of the cluster p{r). Possible bias in our results 
coming from the type of trial wave function used in the 
simulations is resolved with the use of pure estimators 
which ensure unbiased results. (26| 

In Fig. m the density distributions of TJ, clusters hav- 
ing 40, 80, 120, 180, 240 and 320 atoms are plotted. The 
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density profiles show that the cluster size grows when the 
number of atoms in the cluster increases, and that the 
central densities of the largest clusters are very similar 
to the bulk equilibrium density po = 0.007466(7) A'^. 

The increase of cluster size with increasing number of 
atoms can also be seen from the pair distribution function 
P{r) shown in Fig. [5] for the same clusters. P{r) is 
normalized such that J P{r)r^dr — 1. A significant decay 
of the peak height for the largest clusters, as well as the 
growing probability for larger interparticle distances in 
large clusters, is a clear evidence of the size spreading 
tendency. 

The surface thickness of clusters St can be estimated 
from the density profiles as a difference of radii at which 
the central density pc=p{r = 0) has decreased from 90% 
to 10% of its value. From the plotted density profiles in 
Fig. [4] it is obvious that the density error bars are large 
for small distances. In order to determine the central 
density as precisely as possible we have tried to fit the 
density profile with the function used in Rcf. 27j 



p{r) = 



Po 



(1 + e'3('-'^o))<5 



(11) 



where po, /3, rg and 6 are fitting parameters and r is the 
distance to the centre of mass of the cluster. We find that 
for Ti clusters Eq. (HH) can be employed to model den- 
sity profiles of small clusters, while for greater clusters 
the same model reproduces poorly the calculated density 
profiles at small distances. Since the small distances are 
important for our calculation, we have decided to fit the 
calculated density profiles to a constant function for dis- 
tances up to some value ri. We have varied the value of 
ri from 2 A to 4 A, increasing it with the growing size 
of the cluster. We have considered the constant obtained 
with the fit as a central density value and used it in fur- 
ther estimation of the clusters' surface thickness. The 
results for the surface thickness are reported in Table HI 
We can compare our results with the surface thickness of 
"He clusters Sllll. Using the VMC method Pandhari- 



pande et al. [2J] showed that in "He clusters the surface 
thickness is ~7A for clusters > 112. In the case of 
TJ, clusters, for N > 100, the surface thickness is signif- 
icantly greater than in "He (Table |l|). This is expected 
due to the evidently greater interparticle distances in TJ, 
clusters, which is a direct consequence of the shallow at- 
tractive part of T4,-T4, interaction potential. With the 
density functional approach, Stringari et al. calcu- 
lated the surface thickness of several "He clusters and we 
can compare those results with our results for TJ, clusters 
having 20, 40 and 240 atoms. The reported surface thick- 
ness for clusters "He2o, '*He4o and "He24o are respectively 
8.8 A, 9.0 A and 9.3A. A comparison with (T4,)jv reveals 
larger surface thickness of "He at for N— 20,40 and smaller 
surface thickness for A=240. Also, it can be noticed that 
the surface thickness reported by Stringari et al. is not a 
linear function of the number of atoms. Contrary, in the 



case of TJ, clusters, we observe that the surface thickness 
is almost a linear function of A^, up to N= 320 atoms. 
We have tried to predict the surface thickness of clusters 
having more than 320 atoms by fitting our data with 
the function used in Ref i2^ to predict the width of a 
free surface. However, with the present results we have 
not been able to determine the asymptotic value of the 
surface thickness. Saturation should be probably seen 
with results for clusters having more than A=320 atoms, 
but the DMC calculations are already difficult with 320 
atoms. 



CONCLUSIONS 

General characteristics of large spin-polarized tritium 
Tl clusters have been investigated using the DMC ap- 
proach. The ground-state energies of clusters consisting 
of up to 40 TJ, atoms have been compared with previ- 
ously published results. For clusters having more than 
40 atoms the ground-state energies, as well as the struc- 
ture description, are determined for the first time. This 
prediction relies on the use of a very precise potential of 
interaction between T4,-T4, atoms. The present results for 
the ground-state clusters' energy are also used to extract 
the energy per particle of liquid TJ, at equilibrium den- 
sity using a liquid-drop model. The extrapolation using 
a liquid-drop formula gives a value Ey — —3.66(3) K for 
the energy per particle in the equilibrium bulk system, 
which is in very good agreement with the result from a 
recent DMC calculation of the bulk, cq = -3.656(4) AT 

The radii of clusters are calculated with pure estima- 
tors and those results are used to estimate the interval 
in which the unit radius of the liquid is expected. The 
result for the unit radius from the bulk calculation lies 
in the estimated interval. The latter value of the bulk 
unit radius has been employed to estimate the surface 
tension t of bulk TJ,, t=0.08 KA~^. In addition, the sur- 
face thickness of clusters has been estimated from the 
clusters' density profiles. 

As it is already shown for clusters consisting of up to 
40 T\. atoms [ij , it is concluded that large spin-polarized 
tritium clusters are less bound and are more diluted than 
"He clusters with the same number of atoms, i.e., inter- 
particle distances are significantly greater in the corre- 
sponding TJ, clusters. 
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TABLE I: Energy per particle (in K), radii and surface thick- 
ness (in A) of investigated TJ, clusters. 

Figure captions 

FIG. 1: Comparison of calculated ground-state ener- 
gies of clusters (T4,)jv for < 40 atoms (circles) with 
the results reported by Blume et al. in Ref. H (crosses). 
The error bars of the DMC energies are smaller than the 
size of the symbols. 

FIG. 2: Energy per particle E{N)/N for (T^n clus- 
ters reported in TableU Given abscissa is N on an N~^^^ 
scale. The bulk value obtained in is plotted with a 
dashed hne. 

FIG. 3: Radii of (T|)Ar clusters. The abscissa is N on 
an N^^^ scale. The interpolation function ^ is displayed 
with a solid line and the function (fTO|) with a dashed line. 

FIG. 4: Density profiles for several Tj, clusters. Error- 
bars are large for small distances to the centre of mass 
of the clusters, as indicated in the figure for the cluster 
having 240 atoms, and decrease for larger distances. 

FIG. 5: Pair distribution function for several T^ clus- 
ters. 
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FIG. 1: I. Beslic et al. 
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FIG. 2: I. Beslic et al. 




FIG. 3: I. Beslic et al. 
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